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ABSTRACT 

We explore th e em ission properties of a dissipative pulsar magnetosphere model introduced by 
Kalapotharakos et al. (2014), comparing its high energy light curves and spectra, due to curvature 
radiation, with data collected by the Fermi LAT. The magnetosphere structure is assumed to be 
near the force-free solution. The accelerating electric field, inside the light-cylinder, is assumed to be 
negligible, while outside the light-cylinder it rescales with a finite conductivity (a). In our approach we 
calculate the corresponding high energy emission by integrating the trajectories of test particles that 
originate from the stellar surface, taking into account both the accelerating electric field components 
and the radiation reaction forces. First we explore the parameter space assuming different value sets 
for the stellar magnetic field, stellar period, and conductivity. We show that the general properties 
of the model are in a good agreement with observed emission characteristics of young y-ray pulsars, 
including features of the phase resolved spectra. Second we find model parameters that fit each pulsar 
belonging to a group of eight bright pulsars that have a published phase-resolved spectrum. The a 
values that best describe each of the pulsars in this group show an increase with the spin-down rate 
(E) and a decrease with the pulsar age, expected if pair cascades are providing the magnetospheric 
conductivity. Finally, we explore the limits of our analysis and suggest future directions for improving 
such models. 


1. INTRODUCTION 

Pulsars are collapsed cores of massive stars spinning 
at periods in the range icr 3 - 10 sec. They are 
capable of radiating at almost all spectral wavelengths 
(from radio to y-rays), although only a small subset 
have detected y-ray pulsations. Their emission is highly 
anisotropic, originating from particles accelerated along 
the open magnetic field lines by intense electric fields in¬ 
duced by their magnetic field rotation. When pulsars 
display y-ray emission, most of their luminosity is emit¬ 
ted in ~ 100 MeV y—rays. As such they have been one 
of the major class of objects studied by the Fermi y- 
ray Space Tel escope ( Fermi ), usin g the Large Area Tele¬ 
scope (LAT) (Atwood et al. 12 009 ). The release of the 2 nd 
Fermi Pulsar catalog (2PCJ (|Abdo et al.|2013), presents 
light curves and spectra of 117 y-ray pulsars. Of these, 

41 are radio-loud y-ray pulsars discovered through radio 
timing, 36 are new pulsars discovered through their y- 
ray pulsations alone, and a surprising 40 are radio-loud 
millisecond pulsars. We note that for all the radio-loud 
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Fermi pulsars the radio luminosity is orders of magni¬ 
tude smaller than that of y-rays. Nonetheless, since the 
number of radio photons far exceeds the number of the 
y-ray photons radio detection is easier, while y-ray detec¬ 
tion in efficie nt for the nearer and more energetic objects 
(Abdo et al. 2013). 

Pulsar radiation has been studied over the years with 
a variety of models employing different assumptions 
about the location and geome try of emissio n zones , 
such as the Polar Cap (PC) (Daugherty & Hardin; 


1996), the Outer Ga; 


mani & Yadigaroglu 


1983D, IVluslimov V 


IS3 


*1995)' 

Hard! 


(OG) ( |Cneng et al. '1986) , (|Ro- 
the blot Gap (5GJ (|Arons 


Bai & bpitkovsk; 

densky & Spitkovsky) 


yy . . . 

& Kalapotharakos 2010; Kalapotharakos et al. 2012b 

2014). Many models, such as SG and OG implement 


a rdmgl ([2004 ), Separatrix L ayer 
7 2010) and current sheet (Uz- 
mm Petri| |2Q12a| |Contopoulos 


curvature radiation (as in this paper) to describe the 
GeV radiation, while others implement synchrotron ra¬ 
diation powered by the reconnection in the curr ent sheet 
(|Lyubarskii| 1 996| [Uzdensky & Spitkovsky|2014|) . The in¬ 
creased sensitivity of the Fermi LA1\ with its first obser¬ 
vations of a simple exponential cutoff of the Vela pulsar 
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spectrum (Abdo et al. 2009), ruled out the PC emission 
model which predicted a super-exponential cutoff, due 
to magnetic pair attenuation near the neutron star sur¬ 
face. This result, later confirmed for a number of other 
pulsars, established that the pulsar 7 -ray emission origi¬ 
nates in the outer magnetosphere, with the OG and the 
SG models better fitting observations. However, none of 
these outer magnetosphere emission models can currently 
account for the e ntire pulsar lig ht curve phe nomeno logy 
(for example see |Romani & Watters||2010| [Pierbattista 

. 2012, 


et al. 


2014). Furthermore, none of these models are 


consistent with global magnetosphere properties, such as 
current closure. 


trodynamics (FFE) (Spitkovsky 2006 

pn; — — — — 

Kalapotharakos & 

Contopoulos 2009 Petri 2012b), Magnetohydrodynam- 

ical (Komissarov 

2006; Tchekhovskoy et al. 2013) and 

Particle in Cell ( 

Philippov & Spitkovsky 20141 Cerutti 

et al. 

2014 

Chen & Beloborodov 2(JT4) pulsar magne- 


tospnere models gave an impetus to the study of pulsar 
magnetospheres by providing a more realistic picture of 
their outer field geometry than the vacuum model pre¬ 
viously assumed. However, because the ideal MHD FFE 
regime precludes the existence of electric fields parallel 
to the magnetic field (Ey = 0) and thus the acceleration 
of particles and the emission of radiation, more realistic 
dissip ative M HD magnetosp here mo dels hav e been devel¬ 
oped (Kalapothar^^etMj201^ Li et al.| 2012). These 
models allow Tor Ejj 7 ^ 0 and therefore can accommo¬ 
date the production of radiation. They are constructed 
by numerically evolving the time dependent Maxwell’s 
equations, while at the same time employing a form of 
Ohm’s law to relate the current density J to the fields. 
This has the general form 


J = cp 


E x B 

eTTb* 


<tE|| 


( 1 ) 


where J is the current density, E and B are the elec¬ 
tric and magnetic field, and Eq is defined by Bq — Eq = 
B 2 — E 2 , EqBq = E • B, Eq > 0. Eq is a term that 
prevents the drift velocity and hence the current from 
becoming superluminal, while a is a phenomenological 
conductivity, used to relate Ey to the current density J, 
measured in units of the pulsar rotation frequency £ 2 , the 
fundamental frequency in the problem. The macroscopic 
spatial a distribution should be consistent with the un¬ 
derlying microphysics. However, the dependence of the 
conductivity on the microphysics is still not understood. 
Thus, the goal of the recent study of global dissipative 
MHD magnetosphere models is to explore the main prop¬ 
erties of various a values an d distributions. In our first 
studies (Kalapotharakos et al. 20 12 a|c ) a conductivity 
that is mainly uniform (in space) and consta nt in time 
was assumed. In Kalapotharakos et al. (2014), we have 
started to implement a broad range of conductivity val¬ 
ues along the open magnetic field lines while the FFE 
(in reality highly conductive) condition was applied in 
the ‘closed’ region. 

We have employed these dissipative magneto- 
sp here mo dels to generate model 7 -ray light-curves 
(Kalapotharakos et al. 2012c 2014) due to curvature 
radiation (CR). To this end, we calculate the trajectories 
and Lorentz factors 7 l of radiating particles, under 


the influence of both the accelerating magnetospheric 
electric fields and CR-reaction. This approach allowed 
us to relate the observed pulsar emission properties to 
those of the model magnetospheres. Such observables 
are the separation (in phase) of the 7 —ray light curve 
peaks A and also the lag S between the radio emission 
and the fi rst peak of the 7— ray light curve. It was found 
(Kalapotharakos et al. 2014) that matching the model 
d — A distr ibution with that obtained in the 2PC (Abdo 
eFaL||20 13) cannot be achieved with a a constant across 
the entire magnetosphere but requires essentially FFE 
(infinite conductivity, Ey =0) conditions interior to the 
light cylinder (LC), with radius Rlc = c/fi, and large 

(a ~ 30fl) but finite conductivity outside the LC giving 
these models the nomenclature FIDO (FFE Inside, 
Dissipative Outside). 

In the current paper we expand our study by con¬ 
sidering _t he _FIDO models that were presented in 
([Kalapotharakos et al. 2014) by exploring their energetic 
and spectral properties under certain simple assump¬ 
tions. We attempt to compare the detailed luminosities 
and spectra of the dissipative magnetosphere models to 
those obtained by observations. To this end, rather than 
perform a statistical comparison of the FIDO model with 
the entire Fermi population, we decided to compare the 
phase-averaged and phase-resolved spectra predicted by 
the FIDO model with a few very luminous pulsars that 
have published phase-resolved spectra. This being the 
first comparison of spectra predicted by dissipative mod¬ 
els with data, we chose to focus on individual pulsars 
by which we can better understand the model compar¬ 
ison. This comparison will provide the trends and will 
reveal the limitations of this model and its assumptions 
in describing the entire high energy pulsar phenomenol¬ 
ogy depicted by Fermi. 

In §2 we describe the methods employed in our mod¬ 
eling of the pulsar magnetospheres and the production 
of the phase resolved spectra and in §3 we present our 
results. We conclude in §4 with a discussion of their im¬ 
portance of the results and directions to be followed in 
the future studies. 

2. METHODS 

2 . 1 . The Magnetosphere Structure 

The magnetosphere structure is obtained by numeri¬ 
cally solving the time dependent Maxwell equations 


<9B „ „ 

— = -cV x E 

(2) 

<9E 

—— = cV x B — 47rJ 
at 

(3) 


with _a 3D finite_ differ ence time domain technique 
(Kalapotharakos & Contopoulos 2009). We consider the 
presence of a dipole magnetic moment pi at the center of 
the star and that the star itself is a perfect conductor; 
then the boundary condition on the stellar surface for 
the electric field is the FFE condition for a corotating 
magnetosphere: 


E T (Ft x R,) x B — 0 


(4) 


The closure of the system requires a prescription for the 
current density J in terms of the fields (Eq. [l]), as in 
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Figure 1. Examples of photon light curves, between 100 MeV and 50 GeV, produced by a FIDO magnetosphere with P = 0.1s, 
B = 4 • 10 12 G and a = 3012. From left to right are different f (30°, 45°, 60°, 75°, 90°) while from top to the bottom are different a (15°, 
30°, 45°, 60°, 75°, 90°). The overlapping green dots are the phase resolved E cu t. The y axes are rescaled to the relative maximum of the 
light curve, the x axis is rotation phase, from 0 to 1. The values on the y axis are the values of the phase resolved E cu t in GeV. 


Kalapotharakos et al. (2014). When a is very close to 

<9E 

the FFE case (a —» oo,En —0, —» 0), taking the 

O’ dt 

parallel component of Eq. [3] it is possible to approxi¬ 
mate th e elect ric field para lei to the magnetic field as 
Kalap otharakos et al.| (2014) 


E,| = 


c|(V X B) || | ( FFE ) 


47T<J 


(5) 


where the FFE indication implies the corresponding FFE 
value 0 


As in Kalapotharakos et al. 

(|2014), we start with a 

FFE global magnetosphere moc 

el implementing the FFE 

current density J prescription 

Gruzinov||1999) 


J = 


E x B 


cp 


B 2 


c B VxB-E VxE 
47 t B 2 


(6) 


and then compute Ey only outside the LC, located at 
radius Elc = c/fl, using the approximation in Eq. pi of 
the prescription of Eq. [l] for a = 0.3,1, 3, 5,10, 30, oOD 
(a = 3011 was used for the original FIDO model). We did 
not test higher a values since for these the corresponding 
cutoff energies (E cut ) and total y-ray luminosities (E 7 ) 
are well below the ones observed by Fermi. In our FFE 
simulations we use a stellar radius r* = 0.2 Elc- Our 
computational domains extend up to 5Elc and the spa¬ 
tial resolution (i.e. the grid cell size) is O.OIElc- Each 
simulation has been evolved for 2 stellar rotations. 

Here we would like to note that ideal (of infinite resolu¬ 
tion) FFE simulations would lead to infinitesimally thin 


1 This equation comes from a simplified Maxwell equation in 
which the term dE \\/crdt is neglected. 


current sheets and hence to infinitely large V x B || (FFE) 
and E|| values lying within infinitesimally small volumes. 
However, in dissipative solutions (finite a) the thickness 
of the current sheet is expected to be of the order of c/a. 
In our FFE simulations the resolution is finite (O.OIElc) 
and the corresponding current sheet is resolved within a 
few grid cells. This behavior mimics the finite thick¬ 
ness of the current sheet in dissipative solutions. This 
“artificial” FFE thickness of the current sheet is simi¬ 
lar to that corresponding to “real” dissipative solutions 
with a = 3011. Our base FIDO solution is the one with 
a = 3011 presented in Kalapotharakos et al. 2014. We 
have checked that the results corresponding to FFE so¬ 
lutions with a spatial distribution of Ey based on Eq. 5 
(assuming a = 3011) and those of “true” dissipative 
solutions (using explicitly prescription [l]) with uniform 
a = 3011 provide very similar results. 

We also checked the validity of this approximation by 
simulating a few FIDO magnetospheres for a = 0.3, 5fl 
and 10H and various magnetic moment inclination angles 
with respect to the rotational axis, a (fixing the period, 
P, and B); we found that especially for a < 45° and for 
a < 511 the field structure deviates substantially from the 
FFE field structure on which our approximation is based. 
This is primarily because a has been assumed to be uni¬ 
form outside the light cylinder. However, the important 
region for the particle acceleration in the FIDO model is 
near the eq uatorial c u rrent s heet outside the light cylin¬ 
der (Kalapotharakos et al.| 2014). In the case of small 
a values, a more accurate treatment is needed (e.g. ap¬ 
plication of the smaller a values only near the equatorial 
current sheet). This study is beyond the scope of the cur¬ 
rent paper and will be addressed in a future work (see 
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§4). In the rest of the paper we consider the FFE nu¬ 
merical models employing an E\\ that rescales according 
to Eq. ©• This treatment, though simplistic, allows the 
exploration of the main trends and correlations among 
the various parameters taking into account a broad range 
of properties of the observed phenomenology (e.g. light 
curves, spectra). 


2.2. Simulation of the emission 

To simulate the y-ray emission of the model pulsar 
magnetospheres we neglect the inverse Compton (IC) 
and synchrotron radiation (SR) contribution and con¬ 
sider only curvature radiation CR, (Jackson 2001]). In 
some models, where pa rticles are ener g ized by recon - 
nection in current sheet (Uzdensky & Spit kovsky| |2014), 
SR is produced by high temperature particles up to 100 
GeV. In the type of model we are considering here, where 
particles are energized only by the induced electric field 
parallel to the magnetic field, the particles can emit SR 
only if they acquire some pitch angles. In su ch “gap” 
models (|Harding et al.||2008| |Tang et al.||2008) SR does 
not give photons more energetic than 100 IVleV. IC re¬ 
quires a high-altitude non-thermal X-ray emission com¬ 
ponent that is present only for a few Fermi pulsars such 
as the Crab. To save computational time in calculating 
the emission, a random subset of 1.5 • 10 6 initial parti¬ 
cle orbit positions were selected on the pulsar polar cap. 
Their orbits were integrated under the influence of the 
local parallel electric field E\\ and losses due to curva¬ 
ture radiation to compute their Lorentz factor yz, as a 
function of position, in the non-rotating inertial frame. 
Then we selected randomly 5 • 10 10 particles positions 
along these trajectories between the star surface and 2.5 
times the light cylinder radius Rlc, where we computed 
the spectrum of the locally emitted CR and its direction 
(£,</>), taking into account the time delays and assum¬ 
ing that the CR photons are emitted in the direction of 
particle motion (( is the observer’s inclination relative 
to the pulsar rotation axis and <f> the azimuth). For each 
particle we calculated the number of photons emitted per 
second by CR and stored them in a 3D matrix: 100 bins 
of resolution for the azimuth (between 0° and 360°), 
180 bins of resolution for £ (between 0° and 180°) and 
114 energy bins equally spaced logarithmically between 
0.01 and 50 GeV. 

We calculated the luminosity by rescaling the area of 
the polar cap according to the period of the star and to 
a fixed radius of 10 km (in the simulation of the FIDO 
model the radius was at 0.2Rlc) and the flux of particles 
from the surface with the Go ldreich-Julian (GJ) density 


Pgj ([Goldreich & Julian] |1969|)), assuming the polar cap 
is small enough to consider the GJ charge density con¬ 
stant^] The adopted pQj flux assumes that the multi¬ 
plicity of the accelerated particles in the regions with 
high Eu is small (^ 1). This assumption is reasonable 
given that much higher multiplicity leads to screening 
of the E || over distances short enough to prevent fluxes 
>> Pgj of accelerated particles. This 3D grid allows us 
to calculate sky maps, light curves and phased resolved 
spectra. The sky maps are produced by summing up over 
the photon energy and plotting the resulting intensity in 


2 This condition is broken in millisecond pulsars because their 
polar cap size is a significant fraction of their radius. 


photons/s in (£, (f) coordinates. Summing over energy 
for a particular viewing angle ( of the sky map produces 
a light curve. Phase-averaged spectra can be obtained 
by averaging the photon distribution in energy over <p 
for a chosen Finally, phase-resolved spectra are ob¬ 
tained by collecting the spectra in each (£,</>)—bin. We 
applied this procedure to pulsar magnetospheres with a 
number of inclination angles a between the rotation and 
magnetic axes (<a = 15°, 30°, 45°, 60°, 75°, 90°), differ¬ 
ent rotation periods (0.01s, 0.03s, 0.Is, 0.3s, Is) and dif¬ 
ferent magnetic fields at the stellar surface (5 • 10 11 G, 
10 12 G, 2 • 10 12 G, 4 • 10 12 G, 8 • 10 12 G, 3 • 10 13 G) as¬ 


suming each time all the different conductivity values 
(a = 0.3D, ID, 3D, 5D, 10D, 30D, 60D) mentioned in the 
previous section. Thus, we tested (6 x7x5x6 = 1260) 
different magnetosphere configurations. 

The sky m aps and ligh t curves shown in 
Kalapotharakos et al. (20141) were computed for 


the bolometric luminosity^ while here we computed 
sky maps and light curves in photons/s in the energy 
band 0.1 — 50 GeV in order to compare with observed 
Fermi light curves and spectra. Our light curves in 
Figure 1 show some extra secondary peaks that do not 
appear in the bolometric light curves, due to relative 
spectral differences with phase. We generated sky maps 
and light curves in bolo metric lumino sity a nd the y are 
very similar to those of |Kalapotharakos et al.| (|2014|) 
except for small differences in the relative heights of the 
light curve peaks (the maximum difference is <20%). 
This kind of difference is due to differences in sampling 
and/or binning. In addition, we smoothed the profiles 
(with a gaussian smoothing with O = 3.6°) to reduce 
the numerical noise. 


2.3. Comparison with the Fermi pulsars 

We used the data from the 2PC, gathered with three 
years of observations acquired by the LAT on the Fermi 


satellite, except for the phase-resolved spectra that have 

appeared in different papers ( 

DeCesar ( 

2013), Abdo et al. 

(2010a), Abdo et al. (2010b)). These c 

Afferent works use 


data collected over different amounts of time (generally 
less than the 3 years of data used for 2PC) but they all 
analyze photons in the same energy range (between 100 
MeV and 50 GeV like in the 2PC). 

We study some particularly luminous y-ray pulsars, 
in particular the ones that have a published phase re¬ 
solved spectrum. This set of eight pulsars contains: 
J0007+7303 (CTA1 pulsar), J0534+2200 (Crab pulsar), 
J0633+1746 (Geminga pulsar), J0835-4510 (Vela pul¬ 
sar), J1057-5226, J1709-4429, J1836+5925, J1952+3252. 
We identified a candidate light curve and its associ¬ 
ated phase-resolved spectrum produced by our model 
for a combination of <a, £, a that best describe those 
of each pulsar. We made a first selection using the num¬ 
ber of peaks in the light curve. The sample of model 
light curves were generated at different viewing angles 
( (equivalent of the latitude) between 5° and 90° with 
a step of 5°. These 18 light curves were generated for 
each inclination angle and for a = ID, 5D, 30D. This 
grid was fine enough to resolve the single/double peaked 
attribute of the light curves, identifying the best values 

3 y = exp[—cc 2 /© 2 ] 
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of (a,C). We then used simulated magnetosphere mod¬ 
els for the whole range of a noted above, selecting the 
closest P and B values to those of the real pulsars. For 
each combination of a and ( selected we matched the 
cutoff energy E cut from the phase averaged spectra (an 
example in Figure [2| using the sam e fom^as was use d 
to fit the pulsar spectra of the 2PC (Abdo~ eFah][2 Q13) 


d E = K (P 

dE \Eq 


-r 


exp 


E 


Er. 


(7) 



Figure 2. An example of a phase-averaged spectrum produced by 
the model (blue) and fitted to a power law with exponential cutoff 
(green). The red vertical line is the fit value of E cu t • The x axis is 
photon energy in eV, the y axis is a quantity proportional to the 
amplitude. 


We looked for a candidate (considering the whole range 
of a values) with an E cut that was within a factor of 2 of 
that observed. We then matched the separation between 
the gamma-ray peaks (with a tolerance of 0.08 in phase, 
where 1.0 is the period), the evolution of light curves with 
energy predicted by the model, and the flux, in compar¬ 
ison with that observed in the 2PC. The distances used 
are those reported in the 2 PC. Considering the maximum 
error in the distance measures and the maximum varia¬ 
tion in the flux of the simulated light curves with respect 
to and the assumed multiplicity of 1 in the particle 
flux, we accept a flux within a factor of 10 of the observed 
value. These were the primary features we considered in 
selecting candidate light curves and spectra. The other 
features we considered were the “half-power width” of 
the light curve peaks and the trend of the phase-resolved 
E C ut • We did not give as much importance to the spectral 
index T because it is more strongly influenced by possible 
contributions from SR and IC emission. We note that the 
discrepancies in the model flux that we would have ob¬ 
tained considering the observed spectral index (T) value 
in the phase averaged spectrum were always less than 
15%. 


3. RESULTS 

3.1. General properties of the model 


As was reported in Kalapotharakos et al. (2014) the 
FIDO model has the interesting property of reproducing 
the correlati on between 5 and A presented in the 2 PC 
(Abdo et al. 2013). The simulated light curves in FigurefI] 
look similar to the observed light curves, except for some 
additional smaller peaks. In Figure[5]of the Appendix we 
show the luminosity skymaps of a FIDO magnetosphere 
with P = 0.1s, B = 4-10 12 G and a = 511; we also plotted 
the trends in flux and E cut with <a, ( and a. We notice 
that the flux and the E cut generally increase with ( and 
decrease with cr, for fixed a. They also increase with de¬ 
creasing P and increasing B. The physical cause of these 
trends are that increasing ( more closely approaches the 
more energetic main caustic (as we can see in the Ap¬ 
pendix in Figure [5]). Increasing a decreases the value of 
E 11 , the maximum value of the particle Lorentz factor 7 l 


(as we noted in subsection 2 . 1 ) and hence the value of 
E cut . Increasing B also increases E\\ (fixing cr), while 
decreasing P reduces the LC radius and thus on the one 
hand decreases the radius of curvature of the particle 
trajectories (in absolute units) and on the other hand in¬ 
creases the absolute En values (the closer to the star the 
emitting region is the higher the corresponding field val¬ 
ues are). These effects, in general, increase the 7 l values, 
the corresponding photon energies ( E cut ), and the total 
power emitted. Increasing a results in a decrease of E cut 
and of the flux: this is due to both E\\ bein g smaller and 
the GJ charge density (|Goldreich & Julian|(|1969|)) at the 
surface of the polar-cap decreasing with a. reaching 0 at 
a = 90° 0 

In Figure [ 6 j of the Appendix we show the dependence 
of the luminosity of our models with the pulsar period 
and magnetic field (red dots) as well as the observations 
of 2PC (blue dots). In Figure [7] we show the dependence 
of E cut for a pulsar with P = 0.1 sec and B = 2 x 10 12 G 
as a function of <a, £, a while in Fig. [ 8 ] the dependence of 
flux of the same pulsar on the same parameters assuming 
a distance of 1 kpc. Generally, the model E cut values for 
a > 3D match the observed values, while for a < 3D 
the model E cut values are larger than observed. The 
model flux values for all a are within the observed range, 
except for the largest a values at the smallest a and £. 
Although the model fluxes are about a factor of ten below 
the highest observed fluxes at the largest <a, even for the 
lowest a values, Fig. [ 8 ] is only plotted for one value of P 
and B and does not show model fluxes for the full range 
of these values in the observed population. 

Another feature we investigate is the general trend of 
the phase resolved spectra. From the phase resolved 
spectra that are available it seems that the E cut is al¬ 
ways higher in the 2 nd light curve peak or in the right 
part of a broad single peaked light curves. This behav¬ 
ior is followed by ^ 55% of our model light curves, 23% 
do not show a significant difference in E cut between the 
two peaks, while 22 % exhibit the opposite behavior. So, 
according to the FIDO model, obtaining phase-resolved 
spectra from more pulsars should show some that do not 
follow the trend of higher E cut in the second peak. 


b= 1. 


This expression is the same as used in ( Abdo et al.|2 013} with 


5 For a = 90° we substituted the GJ density at 85° to avoid 
having zero charge density for this case. 
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In a few of the magnetosphere models, the luminos¬ 
ity exceeds E. Since this cannot happen physically, the 
problem is linked to the following cause. The approxima¬ 
tion of Eq.[5lforces the E\\ to change linearly with respect 
to cr. NoneTneless, we found through full simulations us¬ 
ing low a that this process is not exactly linear and the 
E\\ saturates for low a values. Thus, whenever the E\\ 
is higher than it should have been this can contribute to 
an overestimation of the luminosity. Finally, this effect 
is probably due to the overestimation of the distances of 
the given pulsars. 


3.2. Model 7 -ray pulsar identification 

First we performed a blind-search to find the parame¬ 
ters giving light curves and spectra that best match those 
of the Fermi pulsars, disregarding any constraints from 
other wavelengths. We note that on the one hand the 
light-curve form is the most important feature for the 
determination of each pulsar’s a value and the Earth’s 
( value. On the other hand the E cut and L 1 values de¬ 
termine the corresponding a value. We found from one 
to three candidates per pulsar. But of the 8 pulsars for 
which we found model candidates with this blind-search , 
two did not have a go od match with other w avelength 
constraints (taken from Pierbattista et al. (2014)). Look¬ 
ing for candidates again in the zone indicated by these 
constraints we found one possible candidate per pulsar, 
with a less-than-perfect light curve match but still rea¬ 
sonable. This fact shows that our model is not always 
able to identify in an unambiguous way a candidate that 
satisfies all our criteria simultaneously. We report the ge¬ 
ometry and conductivity information of our candidates 
in Table [l] while the fluxes and the E cut are shown in 

Table H 


Table 1 

FIDO Pulsar Candidates - Geometry. 


Name 

Cf] 

#\°] 

~«IT 

TIT 

a[Q\ 

Radio 

J0007+7303 

- 

- 

30 

73 

10 

quiet 

J0534+2200 

- 

61 

75 

63 

30 

caustic 

J0633+1746 

- 

> 60 

45 

87 

3 

quiet 

J0835-4510 

70 

64 

60 

50 

10 

loud 

J1057-5226 

- 

- 

45 

61 

30 

loud 

J1709-4429 

- 

53 

45 

55 

10 

loud 

J1836+5925 

- 

- 

90 

40 

1 

quiet 

J1952+3252 

- 

- 

75 

85 

10 

loud 


a Pierbattista et aT| |2014j) 


Table 2 

FIDO Pulsar Candidates - Energetic^ 


a E cu t are in GeV, fluxes in 10 10 ergs ! cm 2 


Name 

/706s 
^ cut 

po6s 

77106s 

r i 

E C ut 

r 

F-, 

cr[f2] 

J0007+7303 

4.7 

1.4 

4.0 

5.0 

0.9 

0.8 

10 

J0534+2200 

4.2 

1.9 

13 

3.4 

1.3 

16 

30 

J0633+1746 

2.2 

1.2 

42 

3.4 

1.2 

1.5 

3 

J0835-4510 

3.0 

1.5 

91 

4.1 

1.4 

19 

10 

J1057-5226 

1.4 

1.0 

3.0 

1.3 

2.0 

1.5 

30 

J1709-4429 

4.2 

1.6 

14 

4.2 

1.1 

1.1 

10 

J1836+5925 

2.0 

1.2 

6.0 

3.0 

1.1 

0.2 

1 

J1952+3252 

2.5 

1.5 

1.4 

2.2 

1.4 

0.7 

10 


We find a correlation between the cr of the candidates 



Figure 3. The conductivity, cr, for the best candidates as a func¬ 
tion of pulsar characteristic age. The pulsars used in the fit are in 
blue; the one in red is J1057-5226 and was excluded from the fit; 
the solid black line is the power law found by the fit routine (see 
Eq. [9} . 

and their E and the age estimate obtained through the 
period and period derivative (excluding the candidate 
for J1057-5226, since it seems to be an outlier). We fit 
these two correlations with a power law, using the large 
intervals we had on the cr grid as errors: 


y = a + b • x c (8) 

We obtained for the relation a— age: 

a = —1 ± 3, b = 32 ± 23, c = -0.36 ± 0.01 (9) 

While for the relation a — E: 

a = -6 ± 30, b = (7.80 ± 0.02) x 1(T 5 , (10) 

c = 0.146 ±0.004 


both with Xred — 0 - 2 ? where Xred * s the reduced y 2 , that 
highlight the large error bar assumed. These relations 
are shown in Figures [3J [4j These two correlations are 
expected if higher a results from more efficient screen¬ 
ing of E 11 by pairs that are produced in greater numbers 


by younger, more e nergetic pulsars (Daugherty & Hard¬ 


ing 1982 Harding & Muslimov||20lT|. Our 8 candidates 

do not reproduce all of the observed characteristics of 
these pulsars: two had a bit smaller flux than observed 
(within a factor of 15, J1709-4429, and 40, J0633+1746, 
instead of 10 ), two of them do not reproduce the ob¬ 
served bridge emission in the light curve (J0633+1746) or 
the off-peak emission present in the observed light curve 
(J1836+5925). For the candidate of radio quiet pulsar 
Geminga (J0633+1746) having a 0.5 peak separation, 
we interchanged the model peaks in order to match the 
observed relative peak heights, which we are free to do 
since our model light curve has no bridge emission. The 
J0007+7303 candidate has a smaller peak width than ob¬ 
served (the tolerance is ±0.06) and two candidates have 
phase resolved spectra that do not show the observed 
trends (J0007±7303, J1836±5925). In the Appendix, 
we show observed and model light curves and phase- 
resolved spectra for three of the pulsars, Vela, Crab and 
J1952±3252. 
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Figure 4. The conductivity, cr , for the best candidates as a func¬ 
tion of pulsar spin down luminosity E. The pulsars used in the fit 
are in blue; the one in red is J1057-5226 and was excluded from 
the fit; t he s olid black line is the power law found by the fit routine 
(see Eq. |10| >. 

4. CONCLUSION 

We have computed light curves, luminosities, phase- 
resolved and phase-averaged spectra for a large number 
of FIDO models (which were shown to reproduce the 
S — A distribution of 2PC) for different values of a , P, 
P, a and (. All of these emission characteristics of the 
model were compared with light curves, phase-aver aged 
and phase-resolved spectra of eight bright pulsars studied 
with Fermi. In particular we explored whether the com¬ 
mon trends found in the eight published phase-resolved 
spectra, such as a higher E cut in the second peak, are 
absolute predictions of the FIDO model. We found that 
this is in fact a common trend, but not absolute since 
the reverse trend (higher E cut in the first peak) or sim¬ 
ilar E cut in both peaks, is also predicted for a subset of 
parameters. It would therefore be of interest to produce 
phase-resolved spectra for a larger number of pulsars. 

This work shows that the FIDO dissipative model mer¬ 
its further exploration, particularly with respect to the 
correlations seen in Figures [3] and [4[ If these correlations 
are confirmed, it would strengthen the expected connec¬ 
tion between screening of E\\ by efficient pair production 
and high conductivity. 

The major limitations of the present model are the as¬ 
sumption of FFE scaling of E\\ with conductivity outside 
the light cylinder, and the assumption that the magneto¬ 
sphere inside the light cylinder is completely FFE. The 
choice of spatial variation of the conductivity (infinite in¬ 
side and finite outside the LC) was imposed by the need 
to reproduce the pulsar 6-A (radio lag - peak separation) 
correlations. Consistency of the models with these cor¬ 
relations requires no (or very little) emission from within 
the LC, hence the FFE approximation for this region. 
Clearly, emission of radiation requires finite conductiv¬ 
ity and non-zero E\\ somewhere in the magnetosphere, 
necessarily outside the LC. A more realistic distribution 
of a likely requires its variation with distance along the 
field lines and from one line to another. This variation 
prescription should be eventually consistent with the re¬ 


lated microphysical mechanisms (the details of which are 
not well understood, at present at least). Nonetheless, 
the overall distribution cannot be much different than 
that proposed in Kalapotharakos 2014, else it would be 
inconsistent with the y-ray light curve correlations. 

We have also explored the validity of the FIDO model 
by making some test simulations with low a outside the 
light cylinder, still assuming FFE conditions inside. It 
was found that the E\\ in these simulations does not drop 
linearly with a for a < 5£4 and a < 45°, primarily be¬ 
cause there is significant departure from FFE field struc¬ 
ture. Therefore, the fluxes and luminosities are some¬ 
what overestimated for the candidates with smaller a. 
However, we have also begun to explore models with low 
a only near the open field boundary, extending outside 
the light cylinder along the current sheet. This distribu¬ 
tion is more realistic since screening by pair cascades is 
expected for most of the interior field lines. In this case, 
the FFE field field structure is present even for lower a 
values, indicating that the FFE scaling of E\\ is more 
valid for a gap-like distribution of a. 

After having explored the possibilities of FIDO-like 
models, we will continue necessary improvements. First, 
the magnetosphere simulations should be calculated ex¬ 
actly and not assuming the linear behavior of Eq. [5] 
Second, an increase in the grid of simulated values is war¬ 
ranted, in particular for the a values, because it is the 
primary source of error in the determination of the simu¬ 
lated flux. Producing a large atlas with more light curves 
would better track rapid changes through the grid values. 
A further step would be to try to understand the origin 
of the <j parameter and its spatial variation, in particular 
the underlying pair produ ction p rocess at different sites, 


such as the polar caps ( 

Timokhin & Arons|[2013), outer 

gaps ( 

Cheng et al.|l986 

or the current sheet. 
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Figure 5. Skymaps produced by a FIDO magnetosphere with P = 0.1s, B = 4 • 10 12 G and a = 5f2 for the energy range 0.01 — 50 GeV. 
Top panels, left to right, are a = 15°, 30°, middle panels are a = 45°, 60°, and bottom panels are a = 75°, 90°. The scale is normalized 
linearly to 99 percentile of the relative maximum. The x axis is rotation phase </>, in degrees; the y axis is observer angle f, in degrees, with 
respect to the rotation axis. 


We show some figures discussed above illustrating general model properties and the data/model comparison. For 
the latter, the first figure for each pulsar shows the observed and model profile and the superimposed phase resolved 
spectra for T and E cut . The second figure shows the evolution of the observed and model light curve in the energy bins: 
100 MeV 4- 300 MeV, 0.3 GeV 4- 1 GeV, 1 GeV 4- 3 GeV, > 3 GeV, > 10 GeV. On the right are the FIDO predictions, 
on the left the data. 
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Figure 6. FIDO model luminosity over the set of para meters (a. a. P. B ) (red) and the luminosity of the pulsars in the 2PC (blue), 
obtained by assuming the flux correction factor, fa = 1 (Abdo et al. 2013}. The horizontal axes are the logio of periods and the log io 
of magnetic field values at the pole: the first are measured in s , the second in G. The vertical axis is the logio of luminosity measured in 
erg/s. 

VELA - J0835-4510 

The Vela pulsar has P = 0.0895 and B = 3.4 • 10 12 G. The A = 0.43, the peaks “half-height width” is 0.003 and 
0.005. The spectrum has E cut = 3.0 GeV, T = 1.5 and the flux observed is 9.1 • 10 -9 erg/(s • cm 2 ). 

The model candidate has period P = 0.1s and B m 4 • 10 12 G. It is for a = 60°, ( = 50°, cr = 1011. The A = 0.42, 

the peaks “half-height width” is 0.06 and 0.07 where the accuracy is ±0.03 due to the smoothing. Our candidate has 
E cut = 4.1 GeV, T = 1.4 and a flux of 1.9 • 10 -9 erg/(s • cm 2 ). 

CRAB - J0534+2200 

The Crab pulsar has P = 0.0335 and B = 3.8 • 10 12 G. The A = 0.40, the peaks “half-height width” is 0.04 and 0.10. 
The spectrum has E cut = 4.2 GeV, T = 1.9 and the flux observed is 1.3 • 10 -9 erg/(s • cm 2 ). 

The model candidate has period P = 0.035 and B = 4 • 10 12 G. It is for a = 75°, ( = 63°, a = 3011. The A = 0.41, the 

peaks “half-height width” is 0.04 for the first peak where the accuracy is ±0.03 due to the smoothing. The candidate’s 
second peak is too broad or too low. Our candidate has E cut = 3.4 GeV, T = 1.3 and a flux of 1.6 • 10 -9 erg/(s • cm 2 ). 

J1952+3252 

J1952 has P = 0.045 and B = 0.5 • 10 12 G. The A = 0.48, the peaks “half-height width” is 0.04 and 0.1. The 
spectrum has E cut = 2.5GeV, T = 1.5 and the flux observed is 1.4 • 10 -1 ° erg/(s • cm 2 ). 

The model candidate has period P = 0.035 and B = 0.5 • 10 12 G. It is for a = 75°, ( = 85°, a = 10H. The A = 0.5, 
the peaks “half-height width” is 0.06 and 0.06 where the accuracy is ±0.03 due to the smoothing. Our candidate has 
E C ut = 2.2 GeV, T = 1.4 and a flux of 7.4 • 10 -11 erg/(s • cm 2 ). 
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< c 

Figure 7. The evolution of E cu t for a pulsar with P = 0.1s and B = 2 • 10 12 G for different (a, £, cr). On the left there is a linear scale, 
on the right logarithmic. The six panels are the different a starting from the top to the bottom: 15°, 30°, 45°, 60°, 75°, 90°. Line colors 
denote different a (in units of 12), according to the legend at the side of the plot. The x axis is f from 0° to 90°, the y axis is the E cu t 
in GeV. The dashed lines are the maximum and minimum E cu t measured for the Fermi young 7 -ray pulsars. In the linear scale panel the 
minimum E cu t is not plotted. 
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Figure 8. The evolution of the fluxes for a pulsar with P = 0.1s, B = 2 • 10 12 G and distance lkpc for different (a, £, a). The six panels 
are the different a starting from the top to the bottom: 15°, 30°, 45°, 60°, 75°, 90°. Line colors denote different ct (in units of 12), according 
to the legend at the side of the plot. On the x axis different f from 0° to 90°; the y axis is the flux in ergs -1 cm -2 . The dashed lines are 
the minimum and maximum flux measured by Fermi for young 7 -ray pulsars, rescaled to a distance of lkpc. 
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Figure 9. Phase resolved spectra for Vela. Left panels: data from |DeCesar| ([2013). Right panels: model candidate light curves (blue) 
and model spectral index (top) and E cu t (bottom) (green) as a function oFpulsar phase. The horizontal blue line marks the E cu t of the 
phase-resolved spectrum. 
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Figure 10. Light curve energy evolution for Vela. Left: Observed light curve in various energy bands, as noted, from (Abdo et al.|2013|). 
Right: Model light curves for the same energy ranges as the data. 
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Figure 11. Phase resolved spectra for Crab. Left panels: data from DeCesar] |2013|. Right panels: model candidate light curves (blue) 
and model spectral index (top) and E cu t (bottom) (green) as a function of pulsar phase. 
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Figure 12. Light curve energy evolution for Crab. Left: Observed light curve in various energy bands, as noted, from (Abdo et al. 2013J. 
Right: Model light curves for the same energy ranges as the data. 
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Figure 13. Phase resolved spectra for J1952+3252. Left panels: data fromfAbdo et al. (2010a). Right panels: model candidate light 
curves (blue) and model spectral index (top) and E cu t (bottom) (green) as a function of pulsar phase. 
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Figu re 14. Light curve energy evolution for J1952+3252. Left: Observed light curve in various energy bands, as noted, from |Abdo et al. 
2013]). Right: Model light curves for the same energy ranges as the data. 















































